Indice | Previo: Aplicaciones.AgujerosNegros | Siguiente: Aplicaciones.AgujerosNegros.FenomenologiaBasica
Queremos encontrar la métrica en el exterior de una distribución de masa con un total de masa M, con simetrÃa esférica.
Por razones de simetrÃa usaremos coordenadas esféricas $(x^0,x^1,x^2,x^3)=(ct,r,\theta,\phi)$ que en este contexto se conocen como coordenadas de Schwarzschild.
Por estar en el exterior de la masa, la ecuación de campo se reduce a la ecuación en el vacÃo:
$$ R_{\mu\nu}=0 $$Infinitas métricas $g_{\mu\nu}$ son compatibles con esta condiciones. ¿Qué es entonces lo que determina la métrica? Las simetrÃas y condiciones de frontera:
La métrica debe tener simetrÃa esférica. Esta simetrÃa implica:
Suponemos que la métrica es estática. Esto significa que no depende del tiempo (estacionaria) pero además no cambia si se invierte el tiempo $t\rightarrow -t$.
La métrica es asintóticamente plana. Esto implica que a gran distancia $t\rightarrow \infty$ la métrica se reduce a la métrica de Minkowski en coordenadas esféricas:
$$ \mathrm{d}s^2=\mathrm{d}T^2-\mathrm{d}r^2-\mathrm{d}{\ell}^2 $$ donde $T=ct$:
$$ \mathrm{d}{\ell}=r^2\mathrm{d}\theta^2+r^2\sin^2\theta\mathrm{d}\varphi^2 $$ es el elemento de lÃnea sobre la superficie de una esfera.
La métrica más general que cumple esas condiciones es:
$$ \mathrm{d}s^2=\mathrm{e}^{2A(r)}\mathrm{d}T^2-\mathrm{e}^{2B(r)}\mathrm{d}r^2-\mathrm{d}{\ell}^2 $$El uso de la función exponencial tiene como proposito:
Nótese también que la parte tangencial de la métrica no depende de $r$ para garantizar la simetrÃa esférica.
La condición de que la métrica debe ser asintóticamente plana implica que:
$$ \lim_{r\rightarrow \infty}A(r)=1\\ \lim_{r\rightarrow \infty}B(r)=1 $$Para procecer a obtener información sobre las funciones $A(r)$ y $B(r)$ debemos calcular todas las cantidades tensoriales asociadas a la métrica y reemplazar en la ecuación de campo.
Las componentes de la métrica son:
$$ \begin{array}{rcl} g_{TT} & = & \mathrm{e}^{2A(r)}\\ g_{rr} & = & -\mathrm{e}^{2B(r)}\\ g_{\theta\theta} & = & -r^2\\ g_{\varphi\varphi} & = & -r^2\sin^2\theta\\ \end{array} $$Y las de la inversa:
\begin{eqnarray} \nonumber g^{TT} & = & \mathrm{e}^{-2A(r)}\\ \nonumber g^{rr} & = & -\mathrm{e}^{-2B(r)}\\ \nonumber g^{\theta\theta} & = & -\frac{1}{r^2}\\ \nonumber g^{\varphi\varphi} & = & -\frac{1}{r^2\sin^2\theta}\\ \end{eqnarray}El procedimiento es bastante laborioso y para hacerlo se puede utilizar una herramienta de cálculo simbólico como Wolfram Lab Cloud. En el sitio electrónico del libro se deja al estudiante las rutinas y una guÃa para los comandos para realizar este cálculo.
Descargue el archivo de texto plano con las instrucciones de este enlace.
Dado que la métrica es simétrica, los únicos sÃmbolos de Christoffel distintos de cero serán:
\begin{eqnarray} \nonumber \Gamma^{\lambda}_{\lambda\nu} & = & +\frac{1}{2} g^{\lambda\lambda} g_{\lambda\lambda,\nu} \\ \nonumber \Gamma^{\lambda}_{\mu\mu} & = & -\frac{1}{2} g^{\lambda\lambda} g_{\mu\mu,\lambda} \end{eqnarray}Los sÃmbolos de Christoffel diferentes de cero serán:
\begin{eqnarray} \nonumber \Gamma^{T}_{Tr} &= & A^{\prime}\\ \nonumber \Gamma^{r}_{rr} &= & B^{\prime} \\ \nonumber \Gamma^{\theta}_{\theta r} & = & \frac{1}{r} \\ \nonumber \Gamma^{\varphi}_{\varphi r} &= & \frac{1}{r} \\ \nonumber \Gamma^{\varphi}_{\varphi\theta} &=&\cot \theta\\ \nonumber \Gamma^{r}_{TT} &=&A^{\prime} \mathrm{e}^{2(A-B)} \\ \nonumber \Gamma^{r}_{\theta\theta} &=&-r \mathrm{e}^{-2 B} \\ \nonumber \Gamma^{r}_{\varphi\varphi} &=&-\mathrm{e}^{-2 B} r \sin ^{2} \theta \\ \nonumber \Gamma^{\theta}_{\varphi\varphi} &=&-\sin \theta \cos \theta \\ \end{eqnarray}El tensor de Riemann es:
$$ {R^\lambda}_{\mu\nu\kappa}=\{\Gamma^\lambda_{\mu\nu,\kappa} + \Gamma_{\mu\nu}^\delta \Gamma_{\kappa\delta}^\lambda \}-\{\nu\leftrightarrow\kappa\} $$que se puede expandir como:
$$ {R^\lambda}_{\mu\nu\kappa}=\Gamma^\lambda_{\mu\nu,\kappa}-\Gamma^\lambda_{\mu\kappa,\nu}+ \Gamma_{\mu\nu}^\delta \Gamma_{\kappa\delta}^\lambda-\Gamma_{\mu\kappa}^\delta \Gamma_{\nu\delta}^\lambda $$Por las simetrÃa del tensor de Riemann y el hecho que la métrica sea diagonal las únicas componentes distintas de cero son:
$$ {R^{\lambda}}_{\mu\nu\kappa},{R^{\lambda}}_{\mu\nu\lambda} $$es decir:
$$ {R^\lambda}_{\mu\nu\lambda}=\Gamma^\lambda_{\mu\nu,\lambda}-\Gamma^\lambda_{\mu\lambda,\nu}+ \Gamma_{\mu\nu}^\delta \Gamma_{\lambda\delta}^\lambda-\Gamma_{\mu\lambda}^\delta \Gamma_{\nu\delta}^\lambda $$En particular los términos distintos de cero serán:
$$ {R^\lambda}_{\mu\mu\lambda}=\Gamma^\lambda_{\mu\mu,\lambda}-\Gamma^\lambda_{\mu\lambda,\mu}+ \Gamma_{\mu\mu}^\delta \Gamma_{\lambda\delta}^\lambda-\Gamma_{\mu\lambda}^\delta \Gamma_{\mu\delta}^\lambda $$Un laborioso cálculo permite encontrar las componentes diferentes de cero no repetidas:
$$ \begin{array}{l} {R^{T}}_{rrT}=A^{\prime} B^{\prime}-A^{\prime \prime}-\left(A^{\prime}\right)^{2} \\ {R^{T}}_{\theta \theta T}=-r \mathrm{e}^{-2 B} A^{\prime} \\ {R^{T}}_{\varphi \varphi T}=-r \mathrm{e}^{-2 B} A^{\prime} \sin ^{2} \theta \\ {R^{r}}_{\theta \theta r}=r \mathrm{e}^{-2 B} B^{\prime} \\ {R^{r}}_{\varphi \varphi r}=r \mathrm{e}^{-2 B} B^{\prime} \sin ^{2} \theta \\ {R^{\theta}}_{\varphi \varphi \theta }=\left(1-\mathrm{e}^{-2 B}\right) \sin ^{2} \theta \end{array} $$El tensor de Ricci se define como:
$$ R_{\mu\nu}=\sum_\lambda {R^{\lambda}}_{\mu\nu\lambda} $$Dado que los únicos términos distintos de cero del tensor de Riemann son aquellos para los cuales $\mu=\nu$, descubrimos que el tensor de Ricci es diagonal:
$$ R_{\mu\mu}=\sum_\lambda {R^{\lambda}}_{\mu\mu\lambda} $$Para el cálculo por ejemplo de $R_{TT}$, por ejemplo:
$$ R_{TT}={R^r}_{TTr}+{R^\theta}_{TT\theta}+{R^\varphi}_{TT\varphi} $$Para obtener el término ${R^r}_{TTr}$ se debe usar la propiedad:
$$ {R^r}_{TTr}=g^{rr}R_{rTTr}=g^{rr}R_{TrrT}=g^{rr}g_{TT}{R^T}_{rrT} $$Después de un álgebra laboriosa obtenemos:
\begin{eqnarray} \nonumber R_{TT} & = & -\mathrm{e}^{2(A-B)}\left[A^{\prime \prime}+\left(A^{\prime}\right)^{2}-A^{\prime} B^{\prime}+\frac{2 A^{\prime}}{r}\right] \\ \nonumber R_{rr} & = & A^{\prime \prime}+\left(A^{\prime}\right)^{2}-A^{\prime} B^{\prime}-\frac{2 B^{\prime}}{r} \\ \nonumber R_{\theta\theta} & = & \mathrm{e}^{-2 B}\left[1+r\left(A^{\prime}-B^{\prime}\right)\right]-1 \\ \nonumber R_{\varphi\varphi} & = & \sin ^{2} \theta\left\{\mathrm{e}^{2 B}\left[1+r\left(A^{\prime}-B^{\prime}\right)\right]-1\right\} \end{eqnarray}El escalar de curvcatura se calcula usando:
$$ R=g^{\mu\nu}R_{\mu\nu}=g^{00} R_{00}+g^{11} R_{11}+g^{22} R_{22}+g^{33} R_{33} $$Al reemplazar queda:
$$ R=-2 \mathrm{e}^{-2 B}\left[A^{\prime \prime}+\left(A^{\prime}\right)^{2}-A^{\prime} B^{\prime}+\frac{2}{r}\left(A^{\prime}-B^{\prime}\right)+\frac{1}{r^{2}}\right]+\frac{2}{r^{2}} $$Por último el tensor de Einstein es la traza inversa del tensor de Ricci:
$$ G_{\mu\nu}=R_{\mu\nu}-\frac{1}{2}g_{\mu\nu}R $$Dado que ambos $g_{\mu\nu}$ y $R_{\mu\nu}$ son simétricos el tensor de Einstein también lo es:
$$ G_{\mu\mu}=R_{\mu\mu}-\frac{1}{2}g_{\mu\mu}R $$El álgebra en este caso produce:
\begin{eqnarray} \nonumber G_{TT} & = &-\frac{2 \mathrm{e}^{2(A-B)}}{r} B^{\prime}+\frac{\mathrm{e}^{2(A-B)}}{r^{2}}-\frac{\mathrm{e}^{2 A}}{r^{2}} \\ \nonumber G_{rr} & = & -\frac{2 A^{\prime}}{r}+\frac{\mathrm{e}^{2 B}}{r^{2}}-\frac{1}{r^{2}} \\ \nonumber G_{\theta\theta} & = &-r^{2} \mathrm{e}^{-2 B}\left[A^{\prime \prime}+\left(A^{\prime}\right)^{2}+\frac{A^{\prime}-B^{\prime}}{r}-A^{\prime} B^{\prime}\right] \\ \nonumber G_{\varphi\varphi} & = &-r^{2} \mathrm{e}^{-2 B} \sin ^{2} \theta\left[A^{\prime \prime}+\left(A^{\prime}\right)^{2}+\frac{A^{\prime}-B^{\prime}}{r}-A^{\prime} B^{\prime}\right] \end{eqnarray}Con todos los elementos a la mano, podemos ahora preguntarnos por cuál es la condición que debe cumplir $A(r)$ y $B(r)$ para que se cumpla la ecuación de campo:
$$ G_{\mu\nu}=\frac{8\pi G}{c^4}T_{\mu\nu} $$o bien:
$$ R_{\mu\nu}=\frac{8\pi G}{c^4}\left(T_{\mu\nu}-\frac{1}{2}g_{\mu\nu}T\right) $$Como estamos en un punto por fuera del cuerpo, eso implica que ambos:
$$ G_{\mu\nu}=0 $$y
$$ R_{\mu\nu}=0 $$Adicionalmente esta última condición implica también que:
$$ R=0 $$Tenemos un total de 9 condiciones posibles:
\begin{eqnarray} \nonumber -\mathrm{e}^{2(A-B)}\left[A^{\prime \prime}+\left(A^{\prime}\right)^{2}-A^{\prime} B^{\prime}+\frac{2 A^{\prime}}{r}\right] & = & 0\;(R_{TT})\\ \nonumber A^{\prime \prime}+\left(A^{\prime}\right)^{2}-A^{\prime} B^{\prime}-\frac{2 B^{\prime}}{r} & = & 0\;(R_{rr}) \\ \nonumber \mathrm{e}^{-2 B}\left[1+r\left(A^{\prime}-B^{\prime}\right)\right]-1 & = & 0\;(R_{\theta\theta}) \\ \nonumber \sin ^{2} \theta\left\{\mathrm{e}^{2 B}\left[1+r\left(A^{\prime}-B^{\prime}\right)\right]-1\right\} & = & 0\;(R_{\varphi\varphi})\\ \nonumber -2 \mathrm{e}^{-2 B}\left[A^{\prime \prime}+\left(A^{\prime}\right)^{2}-A^{\prime} B^{\prime}+\frac{2}{r}\left(A^{\prime}-B^{\prime}\right)+\frac{1}{r^{2}}\right]+\frac{2}{r^{2}} & = & 0\;(R)\\ \nonumber -\frac{2 \mathrm{e}^{2(A-B)}}{r} B^{\prime}+\frac{\mathrm{e}^{2(A-B)}}{r^{2}}-\frac{\mathrm{e}^{2 A}}{r^{2}} & = & 0\;(G_{TT})\\ \nonumber -\frac{2 A^{\prime}}{r}+\frac{\mathrm{e}^{2 B}}{r^{2}}-\frac{1}{r^{2}} & = & 0\;(G_{rr}) \\ \nonumber -r^{2} \mathrm{e}^{-2 B}\left[A^{\prime \prime}+\left(A^{\prime}\right)^{2}+\frac{A^{\prime}-B^{\prime}}{r}-A^{\prime} B^{\prime}\right] & = & 0\;(G_{\theta\theta}) \\ \nonumber -r^{2} \mathrm{e}^{-2 B} \sin ^{2} \theta\left[A^{\prime \prime}+\left(A^{\prime}\right)^{2}+\frac{A^{\prime}-B^{\prime}}{r}-A^{\prime} B^{\prime}\right] & = & 0\;(G_{\varphi\varphi}) \end{eqnarray}Hay muchas posibilidades pero una de las más sencillas es combinar $\mathrm{e}^{-2A}G_{TT}+\mathrm{e}^{-2B}G_{rr}$ para producir:
$$ \frac{2 \mathrm{e}^{-2 B}}{r}\left(A^{\prime}+B^{\prime}\right)=0 $$que implica obviamente:
$$ A'=-B' $$Está última ecuación se puede integrar para producir:
$$ B(r)=-A(r)+C $$¿Qué papel tiene la constante C? Lo único que hace $C$ es hacer el coeficiente de la componente $r$ de la métrica $\exp[2B(r)]=\exp[-2A(r)+2C]=\exp[-2A(r)]\exp(2C)$ un poco mayor o un poco menor en una constante, sin afectar en nada su dependencia de $r$. Ese factor adicional se puede eliminar con un cambio en las unidades de $r$, tal que en esas unidades $\exp(2C)=1$ y por tanto $C=0$. La relación resultante es entonces:
$$ B(r)=-A(r) $$De modo que la métrica será:
$$ \mathrm{d}s^2=\mathrm{e}^{2A(r)}\mathrm{d}T^2-\mathrm{e}^{-2A(r)}\mathrm{d}r^2-\mathrm{d}{\ell}^2 $$Con esta identidad la ecuación para $G_{TT}$ se convierte en:
\begin{eqnarray} \nonumber -\frac{2 \mathrm{e}^{-4B}}{r} B^{\prime}+\frac{\mathrm{e}^{-4B}}{r^{2}}-\frac{\mathrm{e}^{-2 B}}{r^{2}} & = & 0 \\ \nonumber -2 \mathrm{e}^{-2B}r B^{\prime}+\mathrm{e}^{-2B}-1 & = & 0 \\ \nonumber \frac{\mathrm{d}}{\mathrm{d}r}[r(\mathrm{e}^{-2B}-1)] & = & 0 \\ \end{eqnarray}Esta ecuación se integra trivialmente como:
$$ \mathrm{e}^{-2B}=1+\frac{K}{r} $$donde $K$ es una constante.
De aquà y por $A(r)=-B(r)$ obtenemos:
$$ \mathrm{e}^{2A}=1+\frac{K}{r} $$De modo que arribamos a la métrica más general compatible con las condiciones del problema:
$$ \mathrm{d}s^2=\left(1+\frac{K}{r}\right)\mathrm{d}T^2-\left(1+\frac{K}{r}\right)^{-1}\mathrm{d}r^2-\mathrm{d}{\ell}^2 $$Es posible probar con los programas provistos con la versión electrónica del libro que está métrica tiene tensor de Ricci, tensor de Einstein y escalar de curvatura cero como era de esperarse. Adicionalmente vemos que la métrica es asintóticamente plana, cumpliendo con todas las condiciones impuestas inicialmente.
Ahora bien: ¿cuánto vale la constante $K$?. Para esto nos podemos valer del principio de consistencia.
Sabemos que esta métrica en el lÃmite de campo débil debe ser igual a la métrica e Newton:
$$ \mathrm{d}s^2=\left(1-\frac{R_S}{r}\right)\;\mathrm{d}T^2-\mathrm{d}r^2-\mathrm{d}{\ell}^2 $$donde $R_S=2GM/c^2$ es la escala de longitud asociada a la masa $M$.
Una comparación entre estas dos métricas nos muestra que:
$$ K=-R_S $$y la solución queda finalmente:
$$ \mathrm{d}s^2=\left(1-\frac{R_S}{r}\right)\mathrm{d}T^2-\left(1-\frac{R_S}{r}\right)^{-1}\mathrm{d}r^2-\mathrm{d}{\ell}^2 $$Con la forma final de la métrica obtenida podemos escribir finalmente los sÃmbolos de Christoffel:
\begin{eqnarray} \nonumber \Gamma^{T}_{Tr} &=&\frac{R_S}{2r^2}\frac{1}{1-R_S/r}\\ \nonumber \Gamma^{r}_{rr} &=&-\frac{R_S}{2r^2}\frac{1}{1-R_S/r}\\ \nonumber \Gamma^{\theta}_{\theta r} &=&\frac{1}{r} \\ \nonumber \Gamma^{\varphi}_{\varphi r} &=&\frac{1}{r} \\ \nonumber \Gamma^{\varphi}_{\varphi\theta} &=&\cot \theta\\ \nonumber \Gamma^{r}_{TT} &=&\frac{R_S}{2r^2}\left(1-\frac{R_S}{r}\right)\\ \nonumber \Gamma^{r}_{\theta\theta} &=&-r\left(1-\frac{R_S}{r}\right)\\ \nonumber \Gamma^{r}_{\varphi\varphi} &=&-r\left(1-\frac{R_S}{r}\right)\sin ^{2} \theta \\ \nonumber \Gamma^{\theta}_{\varphi\varphi} &=&-\sin \theta \cos \theta \\ \end{eqnarray}Las componentes del tensor de Riemann son:
\begin{eqnarray} \nonumber {R^{T}}_{rrT} & = & \frac{R_S}{r^3}\frac{1}{1-R_S/r}\\ \nonumber {R^{T}}_{\theta \theta T} = {R^{r}}_{\theta \theta r} & = & -\frac{R_S}{2r} \\ \nonumber {R^{T}}_{\varphi \varphi T} = {R^{r}}_{\varphi \varphi r} = -{R^{\theta}}_{\varphi \varphi \theta } & = & -\frac{R_S}{2r}\sin ^{2} \theta \\ \end{eqnarray}Las ecuaciones geodésicas serán:
\begin{eqnarray} \nonumber \dot{U}_T & = & -\frac{R_S}{r^2}\frac{1}{1-R_S/r}U_T U_r\\ \nonumber \dot{U}_r & = & -\frac{R_S}{2r^2}\left(1-\frac{R_S}{r}\right)U_T^2+\\ \nonumber & & +\frac{R_S}{2r^2}\frac{1}{1-R_S/r}U_r^2 +\\ \nonumber & & +\left(1-\frac{R_S}{r}\right) r U_\theta^2 +\\ \nonumber & & +\left(1-\frac{R_S}{r}\right) r \sin^2\theta U_\varphi^2\\ \nonumber \dot{U}_\theta & = & -\frac{2}{r}U_r U_\theta+\cos\theta\sin\theta\;U_\varphi^2\\ \nonumber \dot{U}_\varphi & = & -\frac{2}{r}U_r U_\varphi-2\cot\theta\;U_\theta U_\varphi\\ \end{eqnarray}donde $U_\mu\equiv\mathrm{d}x^\mu/\mathrm{d}\sigma$.
Para conocer las propiedades de la métrica, podemos programar las ecuaciones de la geodésica e integrarlas numéricamente.
Lo primero es escribir las ecuaciones diferenciales:
def geodesica_schwarzschild(Y,s,Rs):
#Coordenadas y sus velocidades
T,r,q,f,dTds,drds,dqds,dfds=Y
#Factores
A=1-Rs/r
v=Rs/(2*r**2)
#ecuaciones
from numpy import zeros_like,sin,cos
dYds=zeros_like(Y)
#Ecuaciones de velocidades
dYds[:4]=Y[4:]
#Aceleraciones
dYds[4]=-2*v/A*dTds*drds
dYds[5]=-v*A*(dTds)**2+\
v/A*(drds)**2+\
A*r*(dqds)**2+\
A*r*sin(q)**2*(dfds)**2
dYds[6]=-(2/r)*drds*dqds+cos(q)*sin(q)*(dfds)**2
dYds[7]=-(2/r)*drds*dfds-2*cos(q)/sin(q)*(dqds)*(dfds)
return dYds
Definamos las propiedades
# Constantes
c=3e8 # m/s
G=6.67e-11 # m^3 kg^-1 s^-2
M=1.98e30 # kg
# Radio de Schwarzschild
Rs=2*G*M/c**2
# Velocidad caracterÃstica
from numpy import pi,sqrt
r0=10*Rs
v0=sqrt(G*M/r0)
df0ds=(v0/r0/c) # r df/d(ct) = v/c
P=2*pi*c*r0/v0 # Periodo tÃpico
Las condiciones iniciales:
from numpy import pi,array
Y0s=array([
0,r0,pi/2,0.0,
1.0,0.0,0.0,0.9*df0ds
])
from numpy import linspace
ss=linspace(0,3*P,300)
La solución es:
from scipy.integrate import odeint
Ys=odeint(geodesica_schwarzschild,Y0s,ss,args=(Rs,))
Ts=Ys[:,0]
rs=Ys[:,1]
qs=Ys[:,2]
fs=Ys[:,3]
En coordenadas cartesianas:
from numpy import sin,cos
xs=rs*sin(qs)*cos(fs)
ys=rs*sin(qs)*sin(fs)
zs=rs*cos(qs)
Un gráfico de la solución
%matplotlib inline
import matplotlib.pyplot as plt
plt.close("all")
from matplotlib.patches import Circle
fig=plt.figure(figsize=(5,5))
ax=fig.gca()
ax.plot(xs/Rs,ys/Rs,'k-')
cRs=Circle((0,0),1,fill=None,ls='--')
ax.add_patch(cRs)
#Decoración
rang=1.2*max(abs(xs).max(),abs(ys).max())/Rs
ax.set_xlim((-rang,rang))
ax.set_ylim((-rang,rang))
ax.set_xlabel("$x/R_S$")
ax.set_ylabel("$y/R_S$")
ax.xaxis.set_label_coords(0.5,0.05)
ax.yaxis.set_label_coords(0.05,0.5)
ax.grid()
fig.tight_layout()
Para ver una versión animada de esta gráfica vea la versión electrónica de este libro.
Podemos reunir los elementos anteriores usando el siguiente código:
#Condiciones iniciales
r0=20*Rs
#Solución a la ecuación
from numpy import pi,sqrt
v0=sqrt(G*M/r0)
df0ds=(v0/r0/c) # r df/d(ct) = v/c
P=2*pi*c*r0/v0 # Periodo tÃpico
from numpy import pi,array
Y0s=array([
0,r0,pi/2,0.0,
1.0,0.0,0.0,0.8*df0ds
])
from numpy import linspace
Nt=300
ss=linspace(0,3*P,Nt)
#Solución
from scipy.integrate import odeint
Ys=odeint(geodesica_schwarzschild,Y0s,ss,args=(Rs,))
Ts=Ys[:,0]
rs=Ys[:,1]
qs=Ys[:,2]
fs=Ys[:,3]
from numpy import sin,cos
xs=rs*sin(qs)*cos(fs)
ys=rs*sin(qs)*sin(fs)
zs=rs*cos(qs)
#Gráfico
import matplotlib.pyplot as plt
plt.close("all")
from matplotlib.patches import Circle
fig=plt.figure(figsize=(5,5))
ax=fig.gca()
punto,=ax.plot([],[],'bo')
trayectoria,=ax.plot([],[],'k-')
cRs=Circle((0,0),1,fill=None,ls='--')
ax.add_patch(cRs)
#Decoración
rang=1.2*r0/Rs
ax.set_xlim((-rang,rang))
ax.set_ylim((-rang,rang))
ax.set_xlabel("$x/R_S$")
ax.set_ylabel("$y/R_S$")
ax.xaxis.set_label_coords(0.5,0.05)
ax.yaxis.set_label_coords(0.05,0.5)
ax.grid()
fig.tight_layout()
def animacion(i):
punto.set_data(xs[i]/Rs,ys[i]/Rs)
trayectoria.set_data(xs[max(0,i-100):i]/Rs,ys[max(0,i-100):i]/Rs)
return punto,trayectoria
#animacion(0)
from matplotlib import animation
anim=animation.FuncAnimation(fig,animacion,frames=Nt,interval=50,blit=True,repeat=False)
Una versión más controlable de la animación se provee a continuación:
from IPython.display import HTML
from matplotlib import rcParams
rcParams['animation.embed_limit']=2**128
HTML(anim.to_jshtml())
Indice | Previo: Aplicaciones.AgujerosNegros | Siguiente: Aplicaciones.AgujerosNegros.FenomenologiaBasica